set more off

use ops_quantities, clear
keep if year>=2000
sort dist_32_id year
tab year, gen(yr_)


gen prod=q*yield/1000
replace prod=0 if mi(yield )


cap drop arain*
gen arain = l.rain10+l.rain11+l.rain12+rain1+rain2+rain3+rain4+rain5
gen arainsq = arain^2
gen araincu = arain^3

lab var prod "Production"
lab var cas_dc_pos "Casualties"

xtivreg cas_dc_pos i.year  (prod= arain*), fe 
predict e, e

xtivreg2 q (cas_dc_pos =e) yr_*, fe cluster(dist_32_id)

*Nytt

cap drop u
eststo clear

eststo: xtreg cas_dc_pos prod yr_*, fe cluster(dist_32_id)
estadd local estimator OLS
estadd local YFE Yes
estadd local FE Yes
eststo: xtreg prod cas_dc_pos yr_* , fe cluster(dist_32_id)
estadd local estimator OLS
estadd local YFE Yes
estadd local FE Yes

eststo: xtivreg2 cas_dc_pos (prod= l.rain10 l.rain11 l.rain12 rain1 rain2 rain3 rain4 rain5) yr_*, fe cluster(dist_32_id)
estadd local estimator 2SLS
estadd local YFE Yes
estadd local FE Yes
predict u, e 
eststo: xtivreg2 prod (cas_dc_pos=u) yr_* , fe cluster(dist_32_id)
estadd local estimator 2SLS
estadd local YFE Yes
estadd local FE Yes

*Lagged
eststo: xtreg cas_dc_pos l.prod yr_* if year>=2002, fe cluster(dist_32_id)
estadd local estimator OLS
estadd local YFE Yes
estadd local FE Yes
eststo: xtreg prod l.cas_dc_pos yr_* , fe cluster(dist_32_id)
estadd local estimator OLS
estadd local YFE Yes
estadd local FE Yes

eststo: xtivreg2 cas_dc_pos (l.prod= l2.rain10 l2.rain11 l2.rain12 ///
  l.rain1 l.rain2 l.rain3 l.rain4 l.rain5) yr_*, fe cluster(dist_32_id)
estadd local estimator 2SLS
estadd local YFE Yes
estadd local FE Yes
eststo: xtivreg2 prod (l.cas_dc_pos=l.u) yr_* , fe cluster(dist_32_id)
estadd local estimator 2SLS
estadd local YFE Yes
estadd local FE Yes

esttab using tables_figs/ivtab.tex, booktabs replace type label ///
 drop(yr_* o.yr_* _cons) star(* .1 ** .05 *** .01)  ///
 coef(prod "Opium production" cas_dc_pos "Casualties" ///
   L.prod "Opium production, lagged" L.cas_dc_pos "Casualties, lagged") ///
 stats(N r2 estimator FE YFE, labels(N "R$^2$"  Estimator "District FE" "Year FE")) ///
 subst(_ \-) nonotes
